home *** CD-ROM | disk | FTP | other *** search
/ Turnbull China Bikeride / Turnbull China Bikeride - Disc 1.iso / HENSA / MATHS / PLPLOT / PLPLOT.ZIP / src / plmap.c < prev    next >
Encoding:
C/C++ Source or Header  |  1994-08-25  |  8.7 KB  |  302 lines
  1. /* $Id: plmap.c,v 1.2 1994/08/25 04:04:49 mjl Exp $
  2.  * $Log: plmap.c,v $
  3.  * Revision 1.2  1994/08/25  04:04:49  mjl
  4.  * Eliminated unnecessary variable decls.
  5.  *
  6.  * Revision 1.1  1994/07/29  20:26:10  mjl
  7.  * Function to plot map backgrounds, read from the specified map data file.
  8.  * Contributed by Wesley Ebisuzaki.
  9.  *
  10. */
  11.  
  12. /*    plmap.c
  13.  
  14.     Continental Outline and Political Boundary Backgrounds
  15.  
  16.     Some plots need a geographical background such as the global
  17.     surface temperatures or the population density.  The routine
  18.     plmap() will draw one of the following backgrounds: continental
  19.     outlines, political boundaries, the United States, and the United
  20.     States with the continental outlines.  The routine plmeridians()
  21.     will add the latitudes and longitudes to the background.  After
  22.     the background has been drawn, one can use a contour routine or a
  23.     symbol plotter to finish off the plot.
  24.  
  25.     Copyright 1991,1993,1994  Wesley Ebisuzaki
  26.  
  27.     This program can be freely distributed, sale except for media
  28.     costs is prohibited.
  29. */
  30.  
  31. #include "plplotP.h"
  32.  
  33. /*----------------------------------------------------------------------*\
  34.  * void plmap(void (*mapform)(PLINT, PLFLT *, PLFLT *), char *type,
  35.  *            PLFLT minlong, PLFLT maxlong, PLFLT minlat, PLFLT maxlat);
  36.  *
  37.  * plot continental outline in world coordinates
  38.  *
  39.  * v1.4: machine independant version
  40.  * v1.3: replaced plcontinent by plmap, added plmeridians
  41.  * v1.2: 2 arguments:  mapform, type of plot
  42.  * v1.1: change buffersize for amiga - faster reads
  43.  *
  44.  * mapform(PLINT n, PLFLT *x, PLFLT *y) is a routine to transform the
  45.  * coordinate longitudes and latitudes to a plot coordinate system.  By
  46.  * using this transform, we can change from a longitude, latitude
  47.  * coordinate to a polar stereographic project, for example.  Initially,
  48.  * x[0]..[n-1] are the longitudes and y[0]..y[n-1] are the corresponding
  49.  * latitudes.  After the call to mapform(), x[] and y[] should be replaced
  50.  * by the corresponding plot coordinates.  If no transform is desired,
  51.  * mapform can be replaced by NULL.
  52.  * 
  53.  * type is a character string. The value of this parameter determines the
  54.  * type of background. The possible values are,
  55.  *
  56.  *     "globe"        continental outlines
  57.  *     "usa"        USA and state boundaries
  58.  *     "cglobe"    continental outlines and countries
  59.  *     "usaglobe"    USA, state boundaries and continental outlines
  60.  * 
  61.  * minlong, maxlong are the values of the longitude on the left and right
  62.  * side of the plot, respectively. The value of minlong must be less than
  63.  * the values of maxlong, and the values of maxlong-minlong must be less
  64.  * or equal to 360.
  65.  * 
  66.  * minlat, maxlat are the minimum and maximum latitudes to be plotted on
  67.  * the background.  One can always use -90.0 and 90.0 as the boundary
  68.  * outside the plot window will be automatically eliminated.  However, the
  69.  * program will be faster if one can reduce the size of the background
  70.  * plotted.
  71. \*----------------------------------------------------------------------*/
  72.  
  73. #define MAP_FILE ".map"
  74. #define OFFSET (180*100)
  75. #define SCALE 100.0
  76. #define W_BUFSIZ    (32*1024)
  77.  
  78. void
  79. plmap(void (*mapform)(PLINT, PLFLT *, PLFLT *), char *type,
  80.       PLFLT minlong, PLFLT maxlong, PLFLT minlat, PLFLT maxlat)
  81. {
  82.     PLINT wrap, sign;
  83.     int i, j;
  84.     PLFLT bufx[200], bufy[200], x[2], y[2];
  85.     short int test[400];
  86.     register FILE *in;
  87.     char filename[100];
  88.  
  89.     unsigned char n_buff[2], buff[800];
  90.     int n;
  91.     long int t;
  92.  
  93. #ifdef AMIGA
  94.     int error;
  95.     char *io_buffer;
  96. #endif    
  97.  
  98.     /*
  99.      * read map outline 
  100.      */
  101.     strcpy(filename,type);
  102.     strcat(filename,MAP_FILE);
  103.  
  104.     if ((in = plLibOpen(filename)) == NULL)
  105.     return;
  106.  
  107. #ifdef AMIGA
  108.     /* increase the size of the IO buffer = faster reads */
  109.     /* ANSI C -- however haven't found it necessary on UNIX */
  110.     io_buffer = (char *) malloc(W_BUFSIZ);
  111.     if (io_buffer == NULL) {
  112.     fprintf(stderr, "ran out of memory: plmap\n");
  113.     return;
  114.     }
  115.     error = setvbuf(in,io_buffer,_IOFBF,W_BUFSIZ);
  116.     if (error != 0) {
  117.     fprintf(stderr, "plmap: setvbuf err %d\n",error);
  118.     free(io_buffer);
  119.     return;
  120.     }
  121. #endif
  122.  
  123.     for (;;) {
  124.     /* read in # points in segment */
  125.     if (fread(n_buff, sizeof (unsigned char), 2, in) == 0) break;
  126.     n = (n_buff[0] << 8) + n_buff[1];
  127.     if (n == 0) break;
  128.  
  129.     fread(buff, sizeof (unsigned char), 4*n, in);
  130.     if (n == 1) continue;
  131.  
  132.     for (j = i = 0; i < n; i++, j += 2) {
  133.         t = (buff[j] << 8) + buff[j+1];
  134.         bufx[i] = (t - OFFSET) / SCALE;
  135.     }
  136.     for (i = 0; i < n; i++, j += 2) {
  137.         t = (buff[j] << 8) + buff[j+1];
  138.         bufy[i] = (t - OFFSET) / SCALE;
  139.     }
  140.  
  141.     for (i = 0; i < n; i++) {
  142.         while (bufx[i] < minlong) {
  143.         bufx[i] += 360.0;
  144.         }
  145.         while (bufx[i] > maxlong) {
  146.         bufx[i] -= 360.0;
  147.         }
  148.     }
  149.  
  150.     /* remove last 2 points if both outside of domain and won't plot */
  151.  
  152.     while (n > 1) {
  153.         if ((bufx[n-1] < minlong && bufx[n-2] < minlong) ||
  154.         (bufx[n-1] > maxlong && bufx[n-2] > maxlong) ||
  155.         (bufy[n-1] < minlat && bufy[n-2] < minlat) ||
  156.         (bufy[n-1] > maxlat && bufy[n-2] > maxlat)) {
  157.         n--;
  158.         }
  159.         else {
  160.         break;
  161.         }
  162.     }
  163.     if (n <= 1) continue;
  164.  
  165.     wrap = 0;
  166.     /* check for wrap around problems */
  167.     for (i = 0; i < n-1; i++) {
  168.         test[i] = (abs(bufx[i]-bufx[i+1]) > 30.0);
  169.         if (test[i]) wrap = 1;
  170.     }
  171.  
  172.     if (wrap == 0) {    
  173.         if (mapform != NULL) (*mapform)(n,bufx,bufy);
  174.         plline(n,bufx,bufy);
  175.     }
  176.     else {
  177.         for (i = 0; i < n-1; i++) {
  178.         x[0] = bufx[i];
  179.         x[1] = bufx[i+1];
  180.         y[0] = bufy[i];
  181.         y[1] = bufy[i+1];
  182.         if (test[i] == 0) {
  183.             if (mapform != NULL) (*mapform)(2,x,y);
  184.             plline(2,x,y);
  185.         }
  186.         else {
  187.             /* segment goes off the edge */
  188.             sign = (x[1] > x[0]) ? 1 : -1;
  189.             x[1] -= sign * 360.0;
  190.             if (mapform != NULL) (*mapform)(2,x,y);
  191.             plline(2,x,y);
  192.             x[0] = bufx[i];
  193.             x[1] = bufx[i+1];
  194.             y[0] = bufy[i];
  195.             y[1] = bufy[i+1];
  196.             x[0] += sign * 360.0;
  197.             if (mapform != NULL) (*mapform)(2,x,y);
  198.             plline(2,x,y);
  199.         }
  200.         }
  201.     }
  202.     }
  203. #ifdef AMIGA
  204.     free(io_buffer);
  205. #endif
  206. }
  207.  
  208. /*----------------------------------------------------------------------*\
  209.  * void plmeridians(void (*mapform)(PLINT, PLFLT *, PLFLT *), 
  210.  *            PLFLT dlong, PLFLT dlat, PLFLT minlong, PLFLT maxlong, 
  211.  *            PLFLT minlat, PLFLT maxlat);
  212.  *
  213.  * Plot the latitudes and longitudes on the background.  The lines 
  214.  * are plotted in the current color and line style.
  215.  * 
  216.  * mapform(PLINT n, PLFLT *x, PLFLT *y) is a routine to transform the
  217.  * coordinate longitudes and latitudes to a plot coordinate system.  By
  218.  * using this transform, we can change from a longitude, latitude
  219.  * coordinate to a polar stereographic project, for example.  Initially,
  220.  * x[0]..x[n-1] are the longitudes and y[0]..y[n-1] are the corresponding
  221.  * latitudes.  After the call to mapform(), x[] and y[] should be replaced
  222.  * by the corresponding plot coordinates.  If no transform is desired,
  223.  * mapform can be replaced by NULL.
  224.  * 
  225.  * dlat, dlong are the interval in degrees that the latitude and longitude
  226.  * lines are to be plotted. 
  227.  * 
  228.  * minlong, maxlong are the values of the longitude on the left and right
  229.  * side of the plot, respectively. The value of minlong must be less than
  230.  * the values of maxlong, and the values of maxlong-minlong must be less
  231.  * or equal to 360.
  232.  * 
  233.  * minlat, maxlat are the minimum and maximum latitudes to be plotted on
  234.  * the background.  One can always use -90.0 and 90.0 as the boundary
  235.  * outside the plot window will be automatically eliminated.  However, the
  236.  * program will be faster if one can reduce the size of the background
  237.  * plotted.
  238. \*----------------------------------------------------------------------*/
  239.  
  240. #define NSEG 100
  241.  
  242. void 
  243. plmeridians(void (*mapform)(PLINT, PLFLT *, PLFLT *), 
  244.         PLFLT dlong, PLFLT dlat,
  245.         PLFLT minlong, PLFLT maxlong, PLFLT minlat, PLFLT maxlat)
  246. {
  247.     PLFLT yy, xx, temp, x[2], y[2], dx, dy;
  248.  
  249.     if (minlong > maxlong) {
  250.     temp = minlong;
  251.     minlong = maxlong;
  252.     maxlong = temp;
  253.     }
  254.     if (minlat > maxlat) {
  255.     temp = minlat;
  256.     minlat = maxlat;
  257.     maxlat = temp;
  258.     }
  259.     dx = (maxlong - minlong) / NSEG;
  260.     dy = (maxlat - minlat) / NSEG;
  261.  
  262.     /* latitudes */
  263.  
  264.     for (yy = dlat * ceil(minlat/dlat); yy <= maxlat; yy += dlat) {
  265.     if (mapform == NULL) {
  266.         y[0] = y[1] = yy;
  267.         x[0] = minlong;
  268.         x[1] = maxlong;
  269.         plline(2,x,y);
  270.     }
  271.     else {
  272.         for (xx = minlong; xx < maxlong; xx += dx) {
  273.             y[0] = y[1] = yy;
  274.         x[0] = xx;
  275.         x[1] = xx + dx;
  276.          (*mapform)(2,x,y);
  277.         plline(2,x,y);
  278.         }
  279.     }
  280.     }
  281.  
  282.     /* longitudes */
  283.  
  284.     for (xx = dlong * ceil(minlong/dlong); xx <= maxlong; xx += dlong) {
  285.         if (mapform == NULL) {
  286.             x[0] = x[1] = xx;
  287.             y[0] = minlat;
  288.             y[1] = maxlat;
  289.             plline(2,x,y);
  290.         }
  291.         else {
  292.             for (yy = minlat; yy < maxlat; yy += dy) {
  293.                 x[0] = x[1] = xx;
  294.                 y[0] = yy;
  295.                 y[1] = yy + dy;
  296.                 (*mapform)(2,x,y);
  297.                 plline(2,x,y);
  298.             }
  299.         }
  300.     }
  301. }
  302.